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Abstract 

We derive the rate of decay of the tail dependence of the bivariate skew normal 
distribution under the equal-skewness condition 6i = 02;= 0, say. The rate of 
convergence depends on whether 0 > 0 or 0 < 0. The latter case gives rate asymp¬ 
totically identical with the case 0 = 0. The asymptotic behaviour of the quantile 
function for the univariate skew normal is part of the theoretical development. 
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1 Introduction 

The coefficient of lower tail dependence of a random vector X = (Vi, V 2 )^ with marginal 
inverse distribution function and is defined as 

Xl = lim Xl{u), where Xl{u) = P{Xi < F^^{u)\X 2 < F^^ffi)). (1) 

X is said to have asymptotic lower tail dependence if Al exists and is positive. If Az, = 
0, then X is said to be asymptotically independent in the lower tail. This quantity 
provides insight on the tendency for the distribution to generate joint extreme event since 
it measures the strength of dependence (or association) in the lower tails of a bivariate 
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distribution. If the marginal distributions of these random variables are continuous, then 
from ([1]), it follows that \l{u) can be expressed in terms of the copula of X, C{ui, U 2 ), as 


Al(m) = 


PiX,<F,-\u),X2<F,-\u)) 


C{u, u) 


u 


( 2 ) 


F{X2<F^\u)) 

Ramos and Ledford ( 2009 1. continuing the work of Ledford and Tawiil ( 1997 1. studied 


intensively a family of bivariate distributions (which they characterised) which satished 
in particular the condition 


\l{u) = ^L{u) 


(3) 


where L{u) is a slowly varying function (SVF) as u —)■ O’*", and a G (0,1], so that, in fact, 
the value of a could be used for comparison of the degree of tail dependence structure 
between members of the family. The standard bivariate extreme value models correspond 
to a = 1 . _ 

Hua and Joe ( 201l[l developed this idea further and dehned k = 1/a in (E]) as the 
(lower) tail order of a copula. The case 1 < k < 2 is termed as intermediate tail depen¬ 
dence as it represents the copula has some level of positive dependence in the tail but 
not as strong as tail dependence with Al = 0. The tail order n can be used to assess tail 
dependence strength when A^ = 0. 

The standard bivariate normal with correlation coefficient — 1 < p < 1 corresponds 
in ([3]), and hence is an instance of inte rmediate tail dep e ndenc e. A recent 
manifestation of this more or less known result is iu iFuug a.ud Seuetal (120111 1. where it is 
shown that 


to a = 


L{u) ^ 2\ ~*~ ^ (—47r logu) i+p, as u —>■ 0^. 

V 1 -p 

The bivariate skew normal distribution was introduced in Azzalini and Dalla Valle 
(1996) (which is discussed further in Azzalini and Capitanio (1999)). A random vector X is 
said to have a bivariate skew normal distribution distribution, denoted as X ~ 8 X 2 ( 0 , R), 
if the probability density of X is 

/(x) = 202(x, R) $ (6»^x), (4) 


where 0(-,i?) is density of a bivariate normal distribution with mean 0 and correlation 
matrix R and $(■) is the cdf of a univariate standard normal distribution. The correlation 

matrix R and skew vector 0 are dehned as ^ Q, with —1 < p < 1 and 0 = (Oi, O 2 Y' ^ 

respectively. Obviously, the (symmetric) normal is obtained as special case when 0 = 0. 

The results of Lysenko, Roy and Waeber (2009), Bortot (2010) and Padoan (2011) 
show that the skew normal distribution is tail independent, that is: \l = lim„^o+ Al(u) = 

0 . 

The focus of this current note is thus to consider ([3]) in the setting of the skew normal 
distribution. 
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We were motivated to do this not only by interest in generalizing to bivariate skew 
norma l the result for the bivariate normal, but also by the following. iFnng and Seneta 
(1201411 chose to view results such as ([3]) as a rate of convergence result to the limi t value 
Al = 0. The skew t distribution is tail dependent i.e. > 0 and iFung and Senetal (1201411 
obtained for it the rate of convergence result: 


\l{u) — Xl = K,{t], R, 6)u^ + 0{u'i) 


(5) 


as M —)■ 0+, where JC{ri,R,0) is a constant which depends on the distribution’s degrees 
of freedom n, scale matr ix R and skewness vector 6 . This gener alised the results of 
Manner and Seg^ ( 2011 1 and Chicheportiche and Bouchaud ( 2012 1. who had taken the 
rate of convergence standpoint and showed its practical importance, for the (symmetric) 
Student’s t distribution. 

Now, this skew t approaches the skew normal and there are cases where Xl and 
IC{ri,R,0) —)■ 0 as p — >■ cx). But it is clear that the expression for Xl{u) deriving from 
(]5jl, in the limit does not seem to provide useful information about an expression for 
convergence rate for the limit. Our treatment of convergence rate for the skew t case had 
worked in a unihed way when 6 i was not necessarily the same as 62 . For the skew-normal 
we conhned ourselves for this paper to equi-skewness: that is 6*1 = 6*2,= 0. Even so, 
we found that our approach provided qualitatively and quantitavely different results in 
the cases 6^ > 0 and 0 < 0, the latter case being asymptotically (apart from a constant 
multiplier) identical to the symmetric bivariate normal case 6^ = 0. The two cases required 
quite different approaches. 

Upper tail dependence behaviour is expres sed immediately from o ur result on lower 
tail dependence behaviour, using the device in iFung and Senetal (120141) . Section 5. 

The remainder of this paper is set out as follows. In Section 2, we derive the asymptotic 
behaviour of the quantile function for the skew normal which is needed for our subsequent 
proof. In Section 3, we derive the rate of convergence in the form of ([3]) for the skew normal 
distribution. 


2 Asymptotic behaviour of the quantile function 

Under the equal skewness condition i.e. 61 = 62 = 0, both Xi and X 2 will have the same 
marginal distribution. Without loss of generality, we shall focus on Xi. The marginal 
density for Xi is 

fxiixi) = 2 (j){xi)^{Xxi), for xi e R, (6) 

where 0(-) and $(•) are the pdf and cdf of the univariate standard normal and 

A = , + . (7) 

yi + 9^1 - 

This means that both Xi and X 2 have a univariate skew normal distribution with skewness 
parameter A i.e. Xi ~ SN{X), i = 1 , 2 . 
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Using Lemma 2 of Capitanio (2010), we have the inequalities for P{Xi < z), for z < 0, 




TT 


-2 


+ 


_A(1 + A2) 1^a(1 + A2) A3(1 + A2) 


„-4 


< PiXi< 




when A > 0 and 


,-^ 2/2 


TT 


1-1 


-2 


LA(1 + A2)' 


< P{Xi <z) < 


IT |A|(1 + A2) 


Z-\-^"^"l2 - |U-3 


2 




- 22/2 


1-1 


+ ^f|^|-4e-^A2/2 

when A < 0. This means that 


7r|A|(l + A2) 
1 


^-2g-22A2/2 


+ 


|A|(1 + A2)2 |A|3(1 + A2) 


I for A > 0; 

' '2 I -| ^ for A < 0 


( 8 ) 


(9) 


as z ^ —oo. (10) 


Recall, for comparison, the cdf of standard normal has the following asymptotic behaviour 
(see for instance Feller (1968) Chapter VII Lemma 2): 


T / \ 1 II - 

<l)(z) ~ —=\z\ 6 2 ^ as 2 : ^ —00 

V 27r 


( 11 ) 


and neither of the expressions on the right-hand side of fITOl) reduces to (fTTj) by letting 
A ^ 0. 

The quantile results are summarised into the following theorem. The asymptotic 
expressions on the right are given in form convenient for the sequel. 

Theorem 1. Let Xi ~ SN{\), then 


Fi (m) ~ y{u) = 


■\/“ W log(-27rAM log(27rAM)), if\> 0; 
-^-2log |^f^-47rlog(f\/2^)^, z/A<0, 


as M —0^ 


Proof. In order to hnd the asymptotic behaviour of the quantile functions Fi^{u) as 
M —O’*", we shall use Theorem 1 of Fung and Seneta (2011), which requires to hnd —y{u), 
a slowly varying function (SVF) as m ^ 0+, such that Fi{y{u))/u —)■ 1 as m —)■ 0+. For 
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then, according to that theorem, ^(n) ~ y{u) as n —)■ O’*". To proceed, we first note 
that both asymptotic expression on the right of (fTOj) are of the form 




( 12 ) 


where a,c, d > 0, 6 < 0. To solve for g{u) in (IT^ as n ^ O’*" we use the Lambert W 
function. The function itself is defined as the solution w to z = we^ for z > 0: that 
is, as the inverse function of the positive continuous and increasing function xe^, x > 0. 
Moreover, the asymptotic behaviour of Lambert W function is given by 


W(z) = \ogz — log log z + 

logz 

as 2 ; —)■ oo so hL( 2 ;) —>■ cxD as —)■ oo; see Corless et al. (1996). From flT^ . 


(13) 




cd, . £A\a(u)\'^ cd /a\F\b\ 

r|»(“)|‘ ’ =r(«) 


so the LHS of the last expression is in the form of we^. Thus 




(14) 


_ 

where W{-) is the Lambert W function. As u —)■ O’*", u i*’i —>■ cx) and we can combine (11^ 
with (IH]) to get 


cd 

]b\ 




\b\ 


\9iu)\ = 


\b\ \Uy 


log 


- log log 


cd f a\ |6| 
|f)| \u) 


^ a j 1% 


/loglog(w(b" 



cd 


Mwn 

as M —)■ O'*". 

Comparing the expression on the right of flTOl) with flT^ . we have a = 
c = ^(1 + A^) and d = 2 so that ^ 
and d = 2 so that ^ = 1 for A < 0. Substitute these constants into flTB]) to get 


g{u) = -< 


(l + A^) 

log log 


log 


, (15) 

,b = -2, 

for A > 0; and a = b = —1, c = 1/2 

h( 

^ ) 


+0 


(i+Ab ^_ 

2 7v\{l-\-X^)u 

'lQg(^^( ^A(l+A^)n )) 
(1+Ab^_1 


' 7rA(H-A^)Li ^ 


1 + A2 


log(^^( ^A(l+A^)j ) 

I 

log(—27rAM log(27rAM)) 


(16) 
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rsj 


(17) 


1 + A2 

as « —)■ O’*" for A > 0; and 

/ 

g{u) = - <( log 


logM, 


i 1 

u\ TT 


\ 


+ 0 


^loglogij- ^ ^ 




-2 log log(|\/^)^ 

-v/-21ogM, 


(18) 
(19) 

as M —>■ O’*" for A < 0. Now set y{u) as the right-hand side of (fT6|) and (ITSD in the 
respective cases A > 0 and A < 0. It is clear that —y{u) is SVF as n —)■ O'*" from (ITTD and 
ina, and since y{u) —)■ —cx) as n —)• O’*", that Fi{ii{u))/u —)■ 1, using ffT6|) . flTSD and the 
right-hand side of fITUD . □ 

Notice (and compare with f[T5]) i that the asymptotic behaviour of the quantile function 
for the standard normal is 


<h ^(u) ~ 2log(tn\/^i^log^ 

~ 2 logu,, as M —0^. 


( 20 ) 


(See iFnng and Senetal (120111 ) where the above methodology is used to obtain (l20ll .) 


3 Main result 

Similarly to the univariate cdf (and the corresponding quantile function), the rate of 
convergence to zero of the lower and upper tail dependence function depends heavily on 
whether 6^ > 0 or 6^ < 0. The results are summarised into the following theorem. 

Theorem 2. Let X ~ 8 X 2 ( 0 , R) with 61 = 62 = 0. As u ^ O’*", 

(a) if9>0, 

.,3 12 




a 


7rA'^/9(l -|- V TT 


— (27rA) 


l+/3^ 


--- )2[— logu]^ 2 


with A = a = — 


\/i+e2(i-p2) - 


^1+202 (l+p) 


and P = yiV2i+p±S); 


(b) if9<0, 


\ I \ 1+P /iTp. , .-£_ 

Al(m) X —-— w- -(— vrlogn) 

2 v 1 ~ P 


( 21 ) 


( 22 ) 
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Proof. The proof will be divided into two parts depends on whether 6* > 0 or 0 < 0. We 
will first consider the case 6 > 0 and mixture representation forms the basis of this proof. 

For a given pair of {0, S) in (|1]), it has been shown in Azzalini and Capitanio (1999) 
that there is a pair of (a, T) such that X can be represented as a normal mean mixture 
by 


X = cxV + Z, 


(23) 


where 


R6 


OL = 


f g(l+P) 


e(i+p) \ - 


Vl + d'^RO 


and 


\ i/i+202(i+p) ’ y/1+202 (i+p) J {a,a) 


(24) 


T = i? - (1 + ROy^ROO^ R = 


/ i+g2(i_p2-) 

p-e^a-p^) \ 

/ 4 

p-02{l-p2) \ 

1 i+2e^(i+p) 

1+262 (l+p) 1 1 


i+2e2(i+p) 1 

1 P_e2(i_p2) 

i+e2(i_p2) 1 1 

P-e2{i-p2) 

4 

\ 1+202(i+p) 

1+2S2(i+p) y 

\ 1+202(i+p) 

^ / 


Note that 6* > 0 implies that a > 0. This parametrisation also satisfies the condition that 
T is symmetric and positive definite. ~ A^(0, T) is the bivariate normal with 0 mean 
and covariance matrix T, and V ~ Half Normal with pdf 


fviy 



ifv>0; 


= 0, otherwise. 


(25) 


That is H = |hF|, with W ~ A^(0, !)• V is assumed to be distributed independently of 
Z. Obviously, when 0 = 0 (and equivalently ck = 0), we have the usual (symmetric) 
multivariate normal as special case for the distribution of X. 

Continuing, when 6^ > 0 


R(X,<Fy(u),X2<Fy(u)) 

= R(Xi < x,X 2 < x), where x = x{u) = Fy(u) 
= Ev{P{aV + Zi < X, aV + Z2 < x)) 


= Fy{u) 


= Ev\P 


Z\ ^ 2 ; — aV Z 2 ~ 


a/\ a/A ’ a/A a/A 

_ p l pi ^x-aV 

X — aV 


Ev ^max(Zj", Zf) < 


a/X 


where we define 


(26) 


Z* = 


y ~N(0, 


T+ePEEP 


T+epTXpj 

1 
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Using the results from Roberts (1966) and Loperfido (2002), we know that = max(ZJ', Z^) ~ 
SN{f3) i.e. a univariate skew normal distribution with skewness 


/9 = 


^ i+e2(i-p2) _ /(1 — p)(l + 26*2(1 + p)) 


\ 1 + 


i+e2(i-p2) 


(1 + P) 


(27) 


If we combine this fact with the Capitanio bounds in ([H]), (|26ll becomes 

POO 1 , , kO 1 \ 

/ lp-5(l+/3U(2^) 


X 


TT 


X — av 


'0 


TT 


/?(l+/32) 


+ 


fv{v) dv 


/3{l + /3^) /33(l + /?2)y ^ a/X 


X — av 


-4 


fv{v)dv 


< I -P I ^ 


< 


X — aV' 

1 ^_ 1 (i+^ 2 )(™) 2 ^^^ / x-av 
0 tt /3(1 + /32) V a/X 


-2 


fv{v)dv, 


(28) 


which suggests that we need to compare the upper bound in fl28D with 




and we will consider 

r ^ ww) (j^ydvMdv 


= \x\ 


71 


a 


A2/3(l + /32) (l + «^)2 


/v/(n) dn 


Applying integration by parts: 

a^\x\ 1 

a 

7rA2/?(l + /?2) [ 

_A2(1 + 

/■°° a 

2d 

7o A2(l + /32)a;^ 

CK^ [ 

afviO) 


-ew(i+/32)- ^ 


e Pd+/3^)( 2 ^ fv{v) 


X 


dn 


+ 


a 


+ "r)^ 

P(l+/32)a;^ 


dn 


-e “ 


X 




dn 


(1 + "r)^ 


dn 


(29) 




































As 


d !e >jy[v} 

dv 


(l + a^)3 lx 


(1 + "r)^ 




(1 + “r)^ 


[-A2(1 + /3>] 


2 ^-A2(1+/32)(1L.) , 

< ^afv{v) + A^(l + (3‘^)vfv{v) + 


(if 




< oo, 


and 4fv{v) = \/~{~v)e we only need 


ve 


dv 


< oo 


for some ^ > 0 to have dominated convergence in (j2^ and the condition is obvionsly true. 
Thus, 

ir X fv{v) dv ^ ^3j^(o) 


lim 

X^ — OO 


I 1-3 -^(l+/32)a;2 


7rA4/i(l + /i2)2’ 


which also implies that 


1 _i 

—e 

TT 


-|(l+/3^)(^)' 


X 


/3(l + /32) V a/A 


X — OLV 


-2 


fv{v)dv 


a'^ 




7rA^/i(l + V vr 

As the hrst term in the lower bound is the same as the upper bound in 
now consider the higher order term in the lower bound in the form 

/o°° X + ^3(iV)) (^) 

I _Jil(l+/32)2,2 

X °e ^ 


(30) 
we will 


= X 


}_^-^{l+li'^)(-2avx+a'^v'^) ^ f 


+ 




TT 





a 


(- 


+ 


V/i(i+/i2) /d^a+zd^); VAV (i + «r)" 


/y(0) 


(if 


7rA6(l + /i2) V/d(l + /d2) /d3(l + /i2) 

after integration by parts and dominated convergence again as x —)■ — cx). Thus 


X — av 


ie-j0+-3^)(^)’ X ( _?_+_i_I I _ 

!r \l3{l + l3^) I3H1+13^) J \ a/\ 


-4 


fv{v) dv 
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^ ; I |_5 _A,^(1+/32)3;2 

= O \ \x\ e ^ ' 


as X —>■ — CX3. By combining 
P{X^<x,X2<x] 


and (l3TD . we have 


a 




7rA‘^/9(l + 13‘^Y V TT 

as X ^ —cx) and by substituting in ffTB]) we have 

P{X,<F^\u),X^<F^\u)) 

~ ^A4/?(l + vf (-TTV 


xe 


-^ (1+/32) (-—^ log(-27rA« log(27rAu))) 


a 


7rA^/3(l + V tt 

/2 

7rA"^/5(l + V TT 


2 a2(i+/ 32) 1 + y2 ^ ^^_^^ 

(27rA) “^(1 + ^^) (- ) 2 [— log m] “^(1 + ^^) 2 + 


(27rA) 


3 r_ 1 ]/3^ 2 V 


[-log m] 


X U 


as 


e\i + pf 


1 + 


e\i + p^) 


e\i + pf {i + 2e\i + p)\ ,2 


^^(i + A^) = 


1 + 202(1 + p ) \^1 + 02 ( 1 _^ 2 ) 

since A = from ((Tj) and a = from fIMll . Finally, from ([T]) 


1 + 202(1 + p)\^ l + 02(l_p2) 

^ / I —I— /.n I I —t- /; I \ 


P{Xr<F{\u)\X,<F^\u)) = 


P{X^<F^\u),X2<F^\u)) 


a 


2 ,, 11 «2 ,1 + A^ 


-(27rA)^+^ 


u 


)2 [— logu]^ 2 X 


( 31 ) 


7rA^/9(l +/S2)2 Y ^ ^ ^ A 2 

which is (j2Tll and part (a) of the proof is now completed. 

Next, we consider the case 0 < 0. In this part of the proof, we proceed by noting 
that from ([2]) that C{u,u) = so that if = z'^L{z),t > 0 where L{z) is 

a slowly varying function as z —?• O’*", then by (applying with suitable transformation to 
regular variation at 0) a result of de Haan (see Seneta (1976), p. 87), we obtain 


C(u,u) 1 dCiz,z) u^Liu) , 

-= - / —-- dz -—, u 0^. 

u u Jq dz r +1 


(32) 


Therefore, it is sufficient for us to hnd a value of r > 0 which satishes = u^L{u), 

for some slowly varying function L(m), as u —)■ 0+, so that fl52]) holds. 
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From ([T]), using L’Hopital’s rule and some well established basic properties of the 
derivative of copulas (see Nelsen (2006), pp.l3, 41), we have 

dC{u, u) 


du 


= P{X2 < F^\u)\X, = F,-\u)) + P(Xi < Pf («)l^2 = F,-\u)) 
= 2F{X2<x\X^ = x), 


(33) 


by letting x = P^ (u) = P 2 ^(u). By using (jl]) and ([6]), the last expression can be written 
as 

fX 1 „„^2 

2F(X2 < x|Xi = x) = 2 


= 2 


/-oo ^27r(l - p2) 

1 <F(6*-^1 - p^z + 6{1 + p)x) 


_ (x 2 -p^ f <F(6*x + 6 x 2 ) , 
e 2 (i-p 2 ) --—^ 


=e 2 " 


e A = 

that is, when u —)■ O’*". Moreover, 


<F(Ax) 


dz, by letting 2 ; = 


X 2 — px 


Since 0 < 0, we have A = 


4>(Ax) ’ Vl - ■ 

< 0, which implies that <F(Ax) —)• 1 as x —>■ —cx). 


a -p 

— CXD < 2 ; < \l -X 

1 + p 


6{1 — p)x < 6 y/l — p^z < 00 , as 6* < 0; 

29x < 6 y /1 — p^z + 6*(1 + p)x < 00 
<F(26'x) < $(6*^1 — p‘^z + 9{1 + p)x) < 1. 


As a result, 
2 


1 — p 

T+P^ 1 


^..,.^{29x) 


dz < 2 F{X 2 < x|Xi = x) < 2 


/I — p 1 9 

T+t^ 1 e-2^ 


« y/^ ^’(Ax 

Since both $(26'x) and <F(Ax) ^ 1 as x —00 we have 

2 P{X 2 < x|Xi = x) ~ 2 [ 


V^^Xx) 


dz. 


^ e-"F"dz = 2 ^{J- —^x) 


y/^ 


l + p 


so 


2P(X, < Fp{n)\X, = Ff ■(«)) ~ 2<S.(J^^Fp(u)), 


since F^^{u) —00 as u —)■ O’*". Since F^^{u) ~ y{u) as u —)■ O'*" from Theorem 1 

and using (iTTll . we can prove 2<F(^^Pf^(M)) ~ 2 ^{^j^y{u)) as n —)■ 0+ by showing 
[F^^{u)]‘^ — [p(n)]^ —)■ 0, as M —)■ O’*". This is in turn equivalent to showing 


,lFp\FPz))]^-[y(FPz))]^ ^ gP+21og(PM^-4Flog(^x^ ^ 


as 2 ; —>• — 00 . 
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Then 


log( \/-47r log( V2^)) 


e'" X ( X (-4rlog(lzl-'e-‘‘'^)), by lUO 

v/TT 


= 1 + 


Att log \z\ 


27rz‘^ 


—)• 1, as z —)■ —oo. 


This implies that 


2P{X2 < Fi\u)\Xi = F^\u)) - 2<h( 


1 -P 
1 + P 


yiu)) 


1 + P / , N —£_ l:z£ 

(—TrlogMj i+Ptii+p 


as M —O’*", by using flTT|) and flT8|) . Finally, by using 


1 -p 

we have 


F{X2<F,-\u)\X,<F,-\u)) 


1 + p /I + p 
1 -p 


_ P 1 —p 

(—7rlog'u)~i^«^, as M —)■ 0^. 


□ 


The theorem shows that when P < 0 there is minimal difference between symmetric 
and skew normal in terms of the intermediate tail dependence as they share the same 
regular varying index which is On the other hand, when P > 0, it has a larger regular 
varying index by a factor of (1 + 2P^(1 + p)) when compared to the normal case and 
therefore skew normal has smaller intermediate tail dependence than the normal in the 
lower tail. 

We shall finish the paper by briefly discussing the corresponding result for the upper 
tail dependence. 

Corollary 1. Let X ~ SN 2 { 0 ,R) with Q\ = Q 2 = Q- In self evident notation and as 
u —)■ 


(a) if 9 < 0 (so —9 = \9\ > Oj, 

\u{u) = P(Xi > F{\u)\X 2 > F:^\u)) 


(1 -m)^ 


a'" 


-(27rA) 




,1 + A2 


)2 [-log(l -m)]' 


7rA"^/9(l + V TT 

w^thX^ and^ ; 

./i+ebi-pb' ,/i+2e2(i+p] ^ V i+i> ’ 


1 

2 


(b) if 9 > t) (so —9 < t)), 


\u{u) ~ (1 - m)i+p X ^^^Y/^^(-7rlog(l - m)) 1+p. 


Proof. If Y = (^ 1 ,^ 2 )^ = —X where X = (Xi,X 2 )^ with continuous marginal distri¬ 
butions, then from Lemma 1 of Fung and Seneta (2014) (with self-evident notation), we 
have A^(m) = Aj(l — m). The proof is completed by noting that when X ~ 3 X 2 ( 6 , R), we 
have Y = —X ~ SN 2 {— 6 , R) from Azzalini and Dalla Valle (1996) and applying Theorem 


2 . 


□ 
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